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ABSTRACT 

A new formulation of MOND as a modified-potential theory of gravity is pro- 
pounded. In effect, the theory dictates that the MOND potential <j) produced by a 
mass distribution p is a solution of the Poisson equation for the modified source density 
p = — (47rG) _1 V • g, where g = v(\g N \/a,o)g N , and g^ is the Newtonian acceleration 
field of p. This makes <j) simply the scalar potential of the algebraic acceleration field g. 
The theory thus involves solving only linear differential equations, with one nonlinear, 
algebraic step. It is derivable from an action, satisfies all the usual conservation laws, 
and gives the correct center-of-mass acceleration to composite bodies. The theory is 
akin in some respects to the nonlinear Poisson formulation of Bekenstein and Milgrom, 
but it is different from it, and is obviously easier to apply. The two theories are shown 
to emerge as natural modifications of a Palatini-type formulation of Newtonian gravity, 
and are members in a larger class of bi-potential theories. 

Subject headings: galaxies: kinematics and dynamics; cosmology: dark matter, theory. 



1. Introduction 

The only known, full-fledged, nonrelativistic (NR) formulation of MOND has been that of 
Bekenstein & Milgrom (1984). This is a modified gravity theory in which the MOND gravitational 
potential, <fi, produced by a mass density p, is gotten from the nonlinear generalization of the 
Poisson equation: 

V • [/i(|V<A|/a o )V0] = 4irGp. (1) 

Here p(x) is a function characterizing the theory. 1 This field equation is derived from an action, 
and enjoys all the standard conservation laws resulting from the usual space symmetries of the 
underlying action. 2 



1 There are also the obvious generalizations to multi-potential theories, whereby the MOND potential is a sum of 
several potentials, each satisfying an equation such as eq.(l). 

2 Blanchet (2007) has given an interpretation to this theory, as resulting from the omnipresence of a gravitationally 
polarizable medium, and Blanchet & Le Tiec (2008, 2009) extended the idea to an appropriate relativistic version. 
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There is also extensive use in the literature of the pristine formulation of MOND (Milgrom 
1983), in which the MOND acceleration, g, is calculated from the Newtonian value, g^, via an 
algebraic relation of the form 

g = v(9N/a )g N , (2) 

where g N = \g N \. For the case of circular motion in an axisymmetric system, appropriate for 
rotation curve analysis, eq.(2) is, in fact, an exact relation in the class of theories dubbed "modified 
inertia" formulation of MOND (Milgrom 1994a). Expression (2) is also exact in the nonlinear 
Poisson formulation for systems of one-dimensional symmetry, such as a spherical galaxy. In this 
case v(y) is related to fi(x) by v(y) = l/fi(x) , where x and y are related by xfi(x) = y. This 
algebraic formulation is very easy to use-hence its attraction as a wieldy tool-and does capture 
the salient MOND effects in many instances. However, it can definitely not be used as a complete 3 
theory, especially for non-test-particle motion: In the first place g resulting from eq.(2) is, generally, 
not derivable from a potential, with all the adverse effects of this (e.g., there is no conserved 
momentum). Also, when applying eq.(2) to describe the center-of-mass motion of a composite 
system based on the accelerations of its constituents it fails completely. The same is true when 
applying it, e.g., to the galactic external-field effect in the solar system (Milgrom 2009a). This had 
left us with only the nonlinear Poisson formulation of MOND as a reliable, complete theory. 

Here I present a new MOND formulation that combines the benefits of the above two formula- 
tions: It constitutes a complete theory derivable from an action, and enjoying the standard conserva- 
tion laws. Yet, the unavoidable nonlinearity-a direct corollary of the basic tenets of MOND-enters 
in an easy to handle, algebraic manner. The resulting field equations then involve only linear 
differential equations. 

This theory came to light as follows: In a recent paper (Milgrom 2009a), concerning the MOND 
external-field effect in the solar system, I applied the nonlinear Poisson formulation to a specific 
problem: a point mass in a constant background field representing the sun in the field of the Galaxy. 
Beside the exact treatment of the problem, I considered a certain approximation for the MOND 
potential of this configuration. It was defined as the solution of the (linear) Poisson equation with 
a source density p = — (47rG) _1 V • g, with g given by eq.(2). This approximation was justified on 
the basis of p having properties similar to those of the density-like field (47rG) _1 A0, where 4> is the 
exact solution of eq.(l). I have now come to realize that this approximation may, in fact, form the 
basis for a complete MOND theory, standing on its own. 

The nonlinear Poisson formulation of MOND has been applied, analytically and numerically, 
to many problems to which it would be interesting to apply the present theory as well. Among 
these: solar system tests (Milgrom 1986a, 2009a, Bekenstein &: Magueijo 2006), forces on non-test- 
mass bodies (Milgrom 1997,2002a, Dai Matsuo &: Starkman 2008), disc stability and bar formation 
(Brada & Milgrom 1999, Tiret & Combes 2008), two-body relaxation (Ciotti & Binney 2004), 



3 'Complete' in the sense of being applicable to an arbitrary problem in the realm for which it is meant; in the 
present case, to any finite, self gravitating system. 
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dynamical friction (Nipoti k al. 2008), escape speed from a galaxy (Famaey, Bruneton, &; Zhao 
2007, Wu et al. 2007), galaxy interactions (Tiret k Combes 2007, Nipoti Londrillo k Ciotti 2007a), 
and collapse (Nipoti, Londrillo k Ciotti 2007b), triaxial models of galactic systems (Wang, Wu, 
k Zhao 2008, Wu k al. 2009), the external-field effect as applied to dwarf spheroidals and warp 
induction (Brada k Milgrom 2000a,b, Angus 2008), structure formation (e.g., Llinares, Knebe, k 
Zhao 2008), and quite a few more. 

Relativistic extensions of the nonlinear Poisson theory where also propounded (see, e.g., 
Sanders 1997, Bekenstein 2004, 2006, Zlosnik Ferreira k Starkman 2006, 2007, Skordis 2009). The 
bi-, or multi-potential nonrelativistic theories described here have also inspired a class of relativistic, 
bi-, or multi- metric MOND theories (Milgrom 2009c,d, see below). 

The new, quasi-linear MOND (QUMOND) theory-in particular its emergence as a modification 
of a Palatini formulation of Newtonian gravity, and its generalizations-is described in section 2. In 
section 3 I derive some of its general properties pertaining to forces on bodies. Section 4 concern 
the deep-MOND limit of the theory. Section 5 deals with systems in a constant external field, 
relevant for the external field effect, and for the question of the center of mass motion of composite 
systems. Section 6 is a discussion. 



2. Formalism 

Consider a gravitating system of density p(r, t) and velocity field v(r, t). Describe the dynamics 
of the system by the action I = f Ldt, with L = f Cd 3 r, and 

C = C K -C P = "^{2V0 • Vr - alQ[{V4> N /a Q ) 2 ]} + P {\v 2 - 0), (3) 

involving the two potentials (f> and 4> N , where C K = pv 2 /2 is the kinetic energy density, and C P 
is the potential energy density. The density p may be viewed as made up of the masses of the 
constituents p(r, t) = J2i "ii5 3 [r — rj(i)], with each of the masses, rrii, negligible with respect to the 
total mass; so each mass can be treated as a test mass in the field of the rest. 

Varying the action over the particle degrees of freedom gives 

h = -V<j>(ri). (4) 

Varying over <\> gives 

Ac/) N = 4irGp, (5) 

and varying over (p N gives 

A0 = V-KlW iV |/a o )V ( A iV ], (6) 

where v(y) = Q'(y 2 )- Equation (4) tells us that the masses move according to the standard 
Newtonian law of inertia in the potential <f>; thus (p is the MOND potential. Equation (5) tells us 
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that <p N solves the Poisson equation with p as source; thus 4> N is the standard Newtonian potential 
(when we impose the standard boundary condition <p N — > at infinity). Equation (6) tells us that 
the MOND potential satisfies the Poisson equation for the density 

as a source, with = —V(j) N . Thus p would be the density that gives the correct MOND potential 
if we interpret the observations in the framework of Newtonian dynamics. Then, p p = p—p would be 
the phantom mass density, introduced in Milgrom (1986b), which will be interpreted by a Newtonist 
as the density of dark matter. The density p p is what was called in Milgrom (2009a) "the surrogate 
mass density" , and was used as an approximation for the phantom density of the nonlinear Poisson 
theory. 

As usual, the symmetries of the action under space and time translations, and rotations, lead 
to conserved momentum, energy, and angular momentum (see also section 3). 

In effect, this theory starts with the acceleration field defined by the algebraic relation eq.(2), 
and remedies its not being derived from a potential by projecting it on the space of gradient vector 
fields. In other words, write 

v(g N /a )g N = -W-VxA, (8) 

which is a unique decomposition if we require that V</> vanish at infinity, then <j> is the MOND 
potential of the theory. 

Restoration of the Newtonian limit for a — >■ requires Q(z) — > z + Q\ for large z, where Q\ is a 
constant. Space-time scale invariance, which is the defining tenet of the deep-MOND limit a —¥ oo 
(Milgrom 2009b), and the standard normalization of a (defined so that the mass-asymptotic- 
rotational-speed-relation is y 4 = MGa ) for which u(y) ~ y -1 / 2 , dictate Q(z) ~ (4/3)z 3 / 4 + Qo 
for z Cl. Since Q is defined up to an immaterial additive constant, we choose it so that Qq = 0. 

In one dimensional cases-e.g., for spherically symmetric systems, -eq. (6) implies relation (2), 
and all three formulations give the same acceleration field with u(y) = l/p(x), where y = xp(x). 
The QUMOND theory will produce unique solutions if and only if yu(y) is monotonic, which we 
assume; and since this function is y for large values of y, it has to be increasing everywhere, so 
yv' + v > 0. 

For an isolated mass distribution of total mass M, bounded in a finite region, we have asymp- 
totically at infinity <p N w -MG/r. Thus p » (4^G)- 1 (MGa ) 1 /2 r - 2 = (l/3)p M (r/i? M )- 2 , where 
p M = 3M/4:TtR 3 m , R m = (MG/o,,) 1 / 2 , are, respectively, the MOND density and the MOND radius 
for the mass M. It thus follows that the MOND potential, relative to its value at a finite point, 
diverges logarithmically at infinity (as in any MOND theory). Note also that the above asymptotic 
behavior of the fields makes L P logarithmically divergent upon space integration. This is the same 
situation as in the case of the nonlinear Poisson formulation. Differences in L P for systems with 
the same total mass are finite, and only such differences will concern us. 
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We can write the MOND potential, relative to some arbitrary origin r = 0, as 

= ~ G j ^(r) (jrbf " R) = h I M ■ ^ N /a Q )g N ](r') (— ^ " p|) , 0) 

which is finite. This is a closed form expression of the MOND field for an arbitrary (bounded) mass 
distribution, with g^ itself being expressed as usual in terms of p. 

2.1. Emergence from a Palatini formulation 

The proposed modification amounts to replacing (V</>) 2 in the standard Newtonian Lagrangian 
density 

c= ~8^G^ 4>)2+p{ r 2 ~ (t)) (10) 

by 2V(p- — a 2 Q[(V(/>/a ) 2 ], where (f> is an auxiliary potential, called here 4> N because it turns out 
to equal the Newtonian potential. This modification can be performed in two steps: First add to £ a 
term (87TG)" 1 (V0— V^) 2 , which is the same as — (87rG) _1 ((/)— <p)(A(f)— Acf>) up to a divergence. This 
does not modify the theory since varying over <f> gives A(p = A<p, which annuls the new addition to 
the action. The free field part of the new Lagrangian density now has — (V</>) 2 + 2V0 • V</> replacing 
(W) 2 . 

The resulting theory after this step may be viewed as a Palatini-type formulation of Newtonian 
gravity. In the standard Palatini formulation of general relativity (GR), the connection (represent- 
ing an auxiliary gravitational acceleration field) is not assumed, a priori, to be the Levi-Civita con- 
nection of the metric. The analogous approach here would be to introduce an auxiliary-acceleration 
degree of freedom, g, and not assume, a priori, that it is the gradient of the potential. Write then 
the Lagrangian density as 

£ = ^(g 2 -2^V.g) + p(iv 2 -0). (11) 

Variation over g gives g = — V</>, and over </>, Vg = —AirGp, yielding Newtonian gravity. 4 The first 
part in £ falls from the gravitational action — ^(ldnG)" 1 J g 1 ^" 2 'g^" 'R of the relativistic Palatini 
formulation of GR. Then, g descends from the independent connection degrees of freedom, from 
which R is constructed, and 4> from the metric. The second part in C falls from the matter 
Lagrangian density, with matter coupled minimally to the metric, and not to the connection. 
In the relativistic Palatini formulation, extremization over the connection (assumed symmetric) 
establishes it as the Levi-Civita connection of the metric (here, g is established as the gradient of the 
gravitational potential), and extremization over the metric (here, 4>) gives the Einstein equation for 
the connection with the standard matter energy-momentum as source (here, the Poisson equation 
for (f)). 



4 The second term in the first part of C can be replaced by +2V0 ■ g from which it differs by a divergence. 



-6- 



MOND is then introduced in the second step, replacing g 2 in £ of eq. (11) by a%Q(g 2 /a%). To 
get our Lagrangian (3), we now do impose, beyond the standard Palatini formalism, that g is a 
gradient of some auxiliary potential (not of <f> itself): g = — V0; so, (ft, not g, is the basic degree of 
freedom beside 4>. If we do not impose this restriction on g we get instead of eqs.(5)(6) a different 
theory 

Klil/ao)g = -W, V-g = -4vrGp, (12) 
which, quite interestingly, is equivalent to the nonlinear Poisson formulation, eq.(l), in the form 
used, e.g., in Milgrom (1986a) [with g = — yu(|V</>|/a )V</>]. 

To recapitulate, the Lagrangian density (11) underlies a Palatini-like formulation of Newtonian 
gravity, whether or not we restrict g, a priori, to be a gradient. However, the MOND extension, even 
with a given Q(z), does depend on whether we make the restriction or not: With the restriction, 
we get our present QUMOND theory. Without it-varying the action over a larger space of trial 
acceleration fields g-we get the nonlinear Poisson theory. This captures at once the affinity and the 
disparity between the two theories, and also explains why they coincide in cases of one-dimensional 
symmetry, where any vector field is a gradient. 



2.2. Generalizations 

Since the Palatini formulation introduces a second degree of freedom, we can further generalize 
the Lagrangian density 5 (11) by taking 

Cp = " 8^fe ■ F(C ' K ' ~ P<f> ' ^ C = ( ^ /0o)2 ' K = ( W flo)2 ' * = 2 ^ • W°o- ( 13 ) 
Extremization over g and over <f> gives, respectively 

^CS = (14) 

V-(7- K V0 + J- ? g) = 47rG / 9. (15) 

Since by eq.(14) g and V0 must be aligned, one of the variables, e.g. k, determines the other two, 
e.g., by solving the equations = J^k = — which follow from eq.(14). Substituting 

the resulting expressions in eq.(15) then gives So even in this general case 

we get the nonlinear Poisson theory eq.(l), and J-(£,k,£) enters the dynamics of the system only 
through the function \x of a single variable, extracted from it. 

If, however, we assume a priori that g is a gradient field g = — V0, the variety of choices of T 
in the above Lagrangian yields a richer family of theories 

V • {T C V4> - F^4>) = 0, V • (J" K V0 - J- ? V0) = 4vrGp. (16) 



5 We require invariance under (f> — ^ 4> + const. Also we limit ourselves here to actions where only g appears, and 
not its derivative, because we want to compare with theories where g is a gradient of some potential, and we want 
to avoid higher derivatives of that potential. 
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It includes our QUMOND theory with F^ = — KC 1 ^ 2 )) F K = 0, F% = — 1, and also the nonlinear 
Poisson theory with F^ = F% = 0, F K = ^(k 1 / 2 ). 

Another interesting subclass is 

jr =O (()+M( K )-M(( + K + (17) 

[note that ( + k + £ = (V0 - V</>) 2 /a 2 ]. With O = and J\f = k we get Newtonian dynamics. For 
A/" = we have F K = F^ = —M.', Fq = O' — M'; so, defining <p* = cf) — (f>, we get from eq.(16) 

V • {O'[(V4>/a ) 2 ]V4>} = -V • {A4'[(V<f /a ) 2 ]V<f } = 4irGp. (18) 

Thus, (f) = (ft + (ft* , with and <ft* satisfying the nonlinear Poisson equation. This covers the NR 
limit of TeVeS for O oc £, so satisfies the linear poisson equation. 

An interesting subfamily is O = aC, (a ^ 0), M = (3 k, for which the Lagrangian density is 

C = ~8^G W ^ )2 + ~ a o M ^<t> ~ V^) 2 /«o]} + P(^v 2 - <j>), (19) 

or, in terms of (f> and <fi* = (j) — <p: 

C = ~8^G {(a + " 2Q ^^ ' ^* + a( ^* )2 " G o- M [(^^)V«o]} + P(^v 2 - 0), (20) 

where I used the freedom to normalize a, (3, and M so that G is the Newton constant. The field 
equations are then 

V>*(|V^|/ao)V^] = 4ttG P , 
A</> = 47rG7T 1 p + /r x V • (M'V^* ) = V • [(1 - a _1 M')V^], (21) 

with 

^ = /3 _«+^_ M/[( ^ 7ao)2] _ (22) 
a 

One has to solve the nonlinear Poisson equation for cf)*, and then the linear version for 4>. In fact, 
for a + (3 / 0, the second eq.(21) can be written as 

A0=— L-A(a<f + <T), (23) 
a + p 

so the MOND potential is simply a linear combination of the Newtonian potential and a solution 
of the first of eq.(21), which is of the type of eq.(l). 

The case a + (3 = (//* = f3) is an interesting special case as it gives our QUMOND theory 
(see below). 

Consider now the Newtonian and MOND limits of the theory underlaid by the Lagrangian 
density (20). It is enough to consider the spherical case so that applying Gauss theorem to the field 
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equations (21) we get algebraic relations between the (minus) radial accelerations g = d(j)/dr > 
0, g* = d(j)* /dr, and the Newtonian acceleration g N = d(j) N /dr > 0: 

„v=.„ s =(i-«-'MV = ae X+'p)M> s - <24) 

In the Newtonian limit g/g N — > 1 so 

M - A<. = ^Zil. (25) 
a + p — 1 

This relation between a, /3, and .M'^ does not require any tuning between different contributions 
to the action. It reflects our taking G to be the phenomenological Newton constant. Starting with 
a coupling G', we would end up with a relation M' x = a(/3 — G'/G)/(a + (3 — G'/G), which is only 
used to define Newton's constant in terms of a, /3, G', and M' x . 

If M is such that = 0, we must have (3 = 1 to get the correct Newtonian limit (or 

rather G = G'/(3, which allows us to normalize the coefficients so that G' = G, (3 = 1), a is not 
constrained. For M.' x infinite we have to have a> + /3 = 1. When M'^ is finite we have neither /3 = 1 
nor a + (3 = 1. 

We can write g as a function of g N alone, distinguishing between two cases. If a + /3 ^ 

9* = v*(\g N \/a )g N , g = g N + - < ^-—^-g N + —^—- u *(\g N \/a )g N , (26) 

a + p a + p 

where ^*(y), as before, is such that if y = xfi*(x) then x = yu*{y). If a + (3 = we have 

g* = (3- X gN, g = g N + p- 2 {M'[{g N /Pa f] - M'Jg*, (27) 

where here, according to eq.(25), M'^ = (3(/3 — 1). From these we can read the requirements for the 
Newtonian and MOND limits. Consider first the case a + (3 ^ 0. In the MOND limit, g N /a — > 
we have to have g/g N — > (flv/flo) -1 ^ 2 diverging. The last term in expression (26) dominates, and 
we must have [a/ (a + j3)\v*{y)y — Y y 1 ^ 2 - 

In considering the Newtonian limit for the a + j3 ^ case, we have to fork again: If a + f3 = 1 
(which, as we saw, applies when M'^ is infinite) the second term in the second of equations 
(26) vanishes, and so f{y)y has to vanish for y — > oo. However, if the first of eqs.(21) is to 
be elliptic, as it must, yv*(y) must be monotonic function. We saw that in the MOND regime 
[a/ (a + j3)}v*{y)y ps y 1 ! 2 is an increasing function, and so it must be increasing everywhere. This 
means that g/g N — 1 must vanish in the Newtonian limit slower than (^jv/oo)" 1 - This, however, is 
in conflict with solar system constraints (Milgrom 1983, 2009a, Sereno & Jetzer 2006); this would 
thus rule out the case a + (3 = l. 6 



This finding is of the same kind arrived at by Zhao and Famaey (2006) in regard to certain versions of TeVeS. 
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When a + j3 ^ 1, v* has to tend to a constant value for large arguments: v*{y) — >■ (a + /3 — 
l)/a, and this it can do with arbitrary speed without violating the ellipticity condition. However, 
comparing the MOND and Newtonian limits in this case we see that if0<a + /3<l, v* has to 
vanish for some finite value of its argument, and this means that n*{x) = l/u*(y) as to blow up 
at a finite argument value. This is undesirable, and excludes the parameter range < a + /3 < 1, 
leaving us with a + (3 > 1 or a + /3 < 0. There may be other constraints on a, (3 to be investigated 
(e.g., positivity conditions). 

Now return to the interesting boundary case a + (3 = 0. There is no matter-of-principlc 
constraint on how fast M! — M'^ can vanish for large arguments; so this theory can be made to 
approach Newtonian dynamics arbitrarily fast. The MOND limit phenomenology dictates M'(z) ~ 
03/2 z -i/4 fo r z < 1. In fact, we can, when a + (3 = 0, absorb (3 in the definition of <p* , so that 
(/)* — > f3~ 1 (j)*, giving an equivalent theory with Ai(z//3 2 ) — (1 — (3~ 1 )z as the new M.. The new M.' 
then vanishes for high arguments. Without loss of generality we can thus put in this case (3 = 1 
and M.' x = 0, with the equations of motion 

A0*=4ttG/9, = V • [(1 + M')V<j>*] = 4ttG P + V • (M'V </>*). (28) 

This theory is equivalent to the QUMOND theory I started with, with Q(z) = z + M(z), and so 
v(y) = 1 + M'(y 2 ). V(f>* is the Newtonian acceleration, and (47rG) _1 V • (M'V<p*) is the density of 
the "phantom matter" representing DM. In the MOND limit z^-Owe have M.'{z) — >■ z -1 / 4 . 

In all the above examples the theory can be cast as two equations (for two potentials) that 
can be solved separately [as for the system (18)], or sequentially [as for system (21)]. For a general 
T this is presumably not possible and solving the resulting theory is rather more challenging. 

Generalizing even further, we could include several auxiliary acceleration fields, gj, in addition 
to (f), and have C P a function of all the scalars g 2 , (V^>) 2 , g« • gj, gi ■ V(j). Without constraining gj 
to be gradients we still get the nonlinear Poisson theory; with the constraints we seem to get a yet 
richer family. All the above theories are equivalent for a spherical system. 

2.3. Extensions 

Here I describe succinctly two interesting extensions of the class of theories described above. 

These bi-potential theories have inspired a class of relativistic formulations of bi-metric MOND 
(BIMOND) theories (Milgrom 2009c,d). The BIMOND theories not only were constructed in 
analogy with the bi-potential theories discussed here, but, in fact, reduce to them in the NR limit. 
The Lagrangian density of the BIMOND theories is constructed after the fashion of the Lagrangian 
density (19): Instead of two potentials we now have two metrics. One, is the MOND metric, 
descending to the MOND potential in the NR limit, and the other, g^, is an auxiliary one giving 
4> in the limit. Matter couples only to the MOND metric in the standard way, echoing the fact 
that here matter couples only to As in many bimetric formulations discussed in the literature, 
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the free Lagrangian densities (V</>) 2 , (V(^) 2 are replaced by the corresponding Ricci scalars of the 
two metrics. The novelty enters in designing the interaction term between the two metrics. The 
difference between the Levi-Civita connections of the two metrics, 

c; 7 = r^ 7 -f; 7 , (29) 

is a tensor, which, moreover, plays the role of gravitational accelerations. Such acceleration-like 
tensors are obviously crucial in the context of MOND, as they permit us to construct dimensionless 
scalars from a~ l Cp 1 and then take functions of these to serve as interpolating functions between 
the GR and the MOND regime. 7 

I have thus considered actions of the form 

1 = ~vLg J [P9 1/2 R + ug 1/2 R-(gg) 1/4 f(K)atM(r/aZ)}d 4 x + I M ^ (30) 

where in the argument of M, T is a scalar quadratic in the tensor C£ . As in our NR case M! 
plays the role of an interpolating function between the MOND and GR regimes. I M is the matter 
action, with matter degrees of freedom represented by ipi, coupling only to g^ v [g and g are minus 
the determinants of the two metrics, and k = (g/g) 1 / 4 ]. I also allow for 'twin' matter described by 
degrees of freedom Xi, which couples only to g^ v . Since the scalar T contains only first derivatives 
of the metrics this action leads to second order field equations. 

The other interesting extension, inspired, in return, by the BIMOND theories, is the possible 
inclusion of 'twin' matter in our NR theories. This can be effected by adding to our Lagrangian 
density (19) a term 7p(^v 2 — <f>) (which would be the NR limit of I M ). We then have two types of 
matter, each accelerated directly, in the standard manner, only by its own potential, but interacting 
indirectly, 'gravitataionally', through the coupling between their potentials, rather unlike standard 
gravity. The field equations (21) now read 

V • [^*(|V<f |/a„)V<f ] = 4vrG(p - /Wp), 

A<p = AnG^p + /r x V • (A4'V(f ) = V • [(1 - a^M')Vf] + iirGja^p. (31) 

I will discuss the interesting implications of such theories with twin matter in a separate paper 
(Milgrom, in preparation), taking in the rest of this paper p = 0. 

3. Forces and virial relations 

I concentrate hereafter on the simpler, QUMOND special case. Most results are easily carried to 
the more general case. Because of the similarities between the theories, many of the properties of the 
nonlinear Poisson formulation derived in the past (see e.g., Bekenstein & Milgrom 1984, Milgrom 
1986a,b, 1994b, 1997, 1998, 2002a) invite derivation of analogous properties of the QUMOND 
formulation. I derive some of these here and in the following sections. 



7 Here c = 1, otherwise we use the MOND scale-length I = c 2 /a in the dimensionless tensors ICp, 
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3.1. Forces 

The force F on the collection of masses in a volume V is defined as the generator of translations, 
in the sense that under a rigid translation of the masses in V alone by a small increment 5r the 
potential energy changes by 

5L P = -Sr • F. (32) 

Such a translation causes a change in the density bp = — 5r ■ Vp, inside V and zero outside. This 
generates corresponding changes in the potentials; but taking into account the stationarity of L P 
under such changes (since the total mass is fixed, the variation in the potentials vanishes at infinity) 
we have SL P = j v 4>5p = —5r • j v 4>V p. Integrating by parts then gives 



= - / d 3 rpV(/). (33) 
Jv 



It is useful to consider the stress tensor of the gravitational field, P. One way to drive it is to 
write L P as a coordinate scalar in curved space, and consider its variation under a change bg iit in 
the metric, to get P from 

5L P = ^J gWtfrPijSg*, (34) 

where summation over repeated indices is understood. One finds from eq.(3), going back to the 
Euclidean case, 

4vrGP = (g • - ^Q)I + u g N ® g" - g ® g N - g N ® g, (35) 

where I is the unit tensor. 8 For solutions of the equations of motion, the divergence of P is found 
to be 

V • P = ps- (36) 
We can thus write the force on V as an integral over its surface S: 

F= f d 3 rpg= [ P da, (37) 
JV is 

or 

F = ^G i {ugN " g)gJV ' da ~ gNg ' da + (g ' gN ~ Y Q)da - (38) 

The torque on the volume V-the generator of rotations in the above sense-is 

T = / d 3 rpr x g = / r x P • da, (39) 
Jv Je 



8 Since P is derived from the action without reference to whether gN is a gradient, this expression for P, and the 
subsequent expressions that use it, remains valid in the nonlinear Poisson formulation if we replace g N by —/i\7<j}, 
and g by — V<j>. This gives 4-kGP = [/x(V<^>) — ^ Q]I — /A7</> <g> V0, which agrees with the expression given in Milgrom 
(2002a). 
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where I integrated by parts making use of the symmetry of P, or 

T = ^|[rx (»g N - g)}g N ■ da - r x g^g ■ da + (g • - ^-Q)r x da. (40) 

In Appendix A I give other expressions for F and T as surface integrals, which might also be useful. 

For an isolated bounded system we get vanishing total force and torque, as can be seen by 
taking the integration surface in eqs.(38)(40) as the sphere at infinity, and using the asymptotic 
behavior of the fields. This is tantamount to the conservation of the total momentum and angular 
momentum in an isolated system. 



3.2. Virial relations 

Multiplying eq.(6) by 4> N and integrating over all space gives an integral relation satisfied by 
solutions of the field equations: 9 

f ■ g - v{g N /a )gl }d 3 r = 0. (41) 



This relation, together with the useful inequality 

zQ'(z) > ^Q(z), for z > 0, (42) 

implies that the free-field energy (87rG) -1 f d 3 r[2V(f) ■ V(j) N — a 2 Q] is always positive for solutions of 
the field equations. 10 To derive inequality (42) look at H(y) = zQ'(z)-Q(z)/2 = y 2 v(y)-Q(y 2 )/2, 
with z = y 2 . We have H(0) = because Q(0) = 0, and H'(y) = y{yv' + v) > for y > 0, from the 
uniqueness condition; so H(y) > for y > 0. 



Define the virial as 

V = 



= - j d 3 rpr ■ g, (43) 



which is useful in several contexts. For example, since the force density on a system can be written 
as f(r) = p(r)g(r), we can write V = — / d 3 rr ■ f. The virial is the generator of space scaling in 
the sense that under r — > (1 + e)r for small e, L P — > L P + eV (see detailed discussion in Milgrom 
1997). Using eq.(36) we write 



V 



J d 3 r(V ■ P) ■ r = - J d 3 rV • (P • r) + J d 3 rV, (44) 



9 With the MOND behavior of v and the asymptotic fields, each of the two terms in the integrand gives rise to an 
integral that diverges logarithmically; the integral of the difference not only converges, but vanishes. 

10 The energy integral actually diverges logarithmically for an isolated system. 
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where the trace 

V = Tr(P) = —[ V & + g • g" - - a 2 Q], (45) 
Write the first term in eq.(44) as a surface integral to get 

AttGV = J d\[2ug 2 N -\lQ]- 

-J(S ■ S N ~ |2)r -da + ur- g N g N ■ da - r • gg^ • da - r • g^g • da, (46) 

where I made use of the integral relation (41). 11,12 In Appendix A I give another useful expression 
for the virial. 

In a theory in which Q(z) vanishes at small z faster than the MOND behavior Q cx z 3 / 4 , such 
as in Newtonian dynamics, the surface integrals vanish at infinity, and what remains of eq.(46) 
forms the basis of a virial relation upon its integration over time. In the case of MOND both 
the volume integral and the surface integrals are finite. Integrating over a sphere at infinity, all 
the vectors appearing in the surface integral in eq.(46) are radial. It follows then that all cancel 
except (a 2 /2) J E Qr da, which contributes (8irG/3)M(MGa ) 1 / 2 , where M is the total mass of the 
system. 13 So, 

V = ^gJ d ' r[2v9 ^ ~ l a » Q] + l M ( M Ga ) 1/2 . (47) 

If we view p as made up of many discrete small point masses, so that each can be considered 
a test particle in the field of the rest, then the force, Fj, on particle i at position r, is — mjV0(rj) 
(this is not true for a non-test particle). We can then write the virial as 14 

V = -^rj-Fj = -^ miYi ■ fj = ^2 m i w2 i - ^ ^2 m i r i " v i- ( 48 ) 

i i i i 

Use eq.(47) and, as usual, average over a long time, dropping the last term in eq.(48) for systems 
that remain bound over long times, we get for such systems the virial relation 

(v 2 ) ee M- 1 (^mjV 2 ) t = \{MGa^ + ^JL^ J d\\W N - \a\Q\) u (49) 
where ()t denotes the long-time average. 



11 Each of the terms in expression (45) for V behaves as r 3 at infinity, but V itself vanishes faster, so it gives a 
finite integral over space. 

12 The above expressions for F, T, and V are invariant to an addition of a constant to Q, as they should be. 

13 This applies for the particular choice of the additive constant in the definition of Q for which Q(0) = 0. 

14 In the fluid approximation, use of the Euler and continuity equations leads directly to V = / d 3 rpv 2 — 
(d/dt) J d 3 rpr ■ v. 
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There is an analogous relation that holds in the nonlinear Poisson formulation (Milgrom 1994b), 
but the integral there involves the MOND potential, for which we have to solve before we can use 
this virial relation. In the present theory the integral involves only the Newtonian acceleration field 
and is readily calculated from the mass distribution-an example of the added amenability of the 
present formulation. 

4. The deep-MOND limit 

In the deep-MOND limit the Lagrangian density becomes 

c = £fe [( ^ r/a ° )2]3/4 - i^G^ ■ ^ r + p( r 2 - ^ (50) 

Multiplying by a G and absorbing a (j> N — > <j) N gives an action with only a Grrii appearing, as 
required by the basic MOND tenets. Seen differently, under space-time scaling (r,t) — > A(r, t), 
with (j) having zero internal dimension so that 0(r) — > <p(r/X), and <p N having dimension —1, so 
4> N (r) \~ 1 (j) N (r/\), L is invariant. The action is thus multiplied by A; so the equations of motion 
of the masses are invariant. This is the defining requirement of the deep-MOND limit of any theory 
(Milgrom 2009b). 

For the present theory there exists another symmetry of the potential Lagrangian, L P , alone: 
Under space scaling r — > Ar, for a constant A, where, again, <p{r) — > ^(r/A), <p N (r) — > \~ 1 <fi N (r/\), 
the potential Lagrangian density C P (r) — > A~ 3 £ P (r/A), and so L P and the potential action are 
invariant. This means that the field equations in the deep-MOND limit of the theory are invariant 
to space scaling: If <j>(r) is a solution for p(r), then 0(r/A) is a solution for A~ 3 p(r/A). Recall 
that the deep-MOND limit of the nonlinear Poisson equation is invariant to the whole group of 
conformal transformations in space, of which scaling is one (Milgrom 1997). 

4.1. Virial relation 

If the system is everywhere deep in the MOND regime, for which a%Q = (4/3)vg'j f with our 
choice of the additive normalization of Q, the volume integral in eq.(49) vanishes, and we have 

(v 2 ) = \{MGa fl\ (51) 

This is an even more useful virial relation than eq.(49), since it relates the 3-D root-mean-square 
velocity, of any system of test particles deep in the MOND regime, to the total mass, irrespective 
of the mass distribution. 

In many instances we want to view the system as made up of several finite, non-test masses, 
ignoring the internal goings on in each of them, and using only the center of mass velocities of these 
masses (such as galaxies in a small group, or a binary). When the accelerations of these masses are 
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still deep in the MOND regime we can write the above relation (following e.g., the argumentation 
in Milgrom 1997) as 

^) = 2 -(MGa i) )^(l-Y J <ir\ (52) 

i 

where qi = rrii/M. This is based on the relation 

J> • F, = - 2 -M(MGarfl\l - J> 3/2 ), (53) 

i ' i 

which applies in this case. All these results are identical to those for the deep-MOND limit of 
the nonlinear Poisson formulation. The general two-body force in the deep-MOND limit follows 
straightforwardly from this relation, and is thus the same as in the nonlinear Poisson formulation 
(see Milgrom 1997). 

Also discussed there is the relation between the above results and the space scale invariance of 
the nonlinear Poisson formulation. This is also relevant to the present theory. I only note here that 
despite the formal invariance of L P under scaling, the fact that L P itself is infinite leaves room for a 
finite change in L P (which can be viewed as the potential energy of the system) under scaling of the 
density distribution p(r) — > p\(r) = A~ 3 p(r/A). So, under such a change in the mass distribution 
we have 

L P [ PX ] = L P [p] + ln{X v ) = L P [p] + -M{MG ai) ) l / 2 ln{\). (54) 

This violation of the symmetry by the transformation properties of L P is analogous to so called 
"anomalies" that appear in scale- or conformally-invariant field theories, with the virial playing the 
role of anomalous dimension. In a system made of finite, point masses at positions r$ we have 15 

L P (\vi) = L P ( ri ) + ^M(MGa )^ 2 (l - q- /2 )ln(X), (55) 

i 

giving rise to relation (53) by taking its A derivative at A = 1. 



5. A system in a constant external field 

Many times one is dealing with the dynamics in and about a small system that is itself subject 
to the gravitational field of a large, mother system. Examples are the field around a star, a globular 
cluster, or a dwarf companion in a field of a mother galaxy. Instead of solving the full MOND 
subsystem-plus-mother-system problem, we may want to treat the subsystem alone, taking the 
presence of the mother system into account through the constant-acceleration boundary condition 
that it dictates in the vicinity of the subsystem. This approximation is good if we can erect around 



Here, we only move the point masses to new positions Ar;, without dilating them intrinsically by a factor A, as 
would be required for eq.(54) to apply. 
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the subsystem a volume V, whose surface is S, such that (i) V is small enough that the acceleration 
due to the mother system alone is nearly constant in V, (ii) V is large enough so that on E the 
perturbation due to the subsystem is negligible compared with the acceleration due to the mother 
system. Assume then that these conditions are satisfied. Solve for the field of the mother system, 
which gives values and go for the Newtonian and MOND accelerations at the position of the 
subsystem. We do not (and need not) specify any relation between these two accelerations. For a 
spherical mother system they are related by the algebraic relation (2), but in general they are not 
so related, and are not even parallel. 

Our approximation-becoming increasingly better as the mass of the subsystem is small com- 
pared with that of the mother system-is to solve the MOND field equations for the density p of the 
subsystem alone, dictating as boundary conditions at infinity V(p — > —go, —¥ ~Sq, instead of 
the conditions appropriate for an isolated system. We write 

V<T = -g^ + V X) V0 = -g + VV, (56) 

and determine ip from eq.(6) with the boundary condition Vip — >■ at infinity, given that \ ls 
the standard Newtonian potential of the subsystem. Note, importantly, that go drops from the 
equation for the internal potential ip; so, this potential does not depend on go, only on g^. Take 
the direction of the latter to define the positive z axis, and its value in units of a : rj N = \^\/a . 

Unlike the isolated case, where M = f pd 3 r is infinite, for the case of a constant background 
field M is finite, and is given by (Milgrom 2009a) 

M = J vg N ■ da = u (l + *> /3)M, (57) 

where uq and z>o are, respectively, the values of v and its logarithmic derivative at 7] N (—1/2 < Oq < 
0), and M is the mass of the subsystem. 



5.1. Asymptotic behavior of the field 

The asymptotic behavior of p is obtained from that of the Newtonian field, and is given by 

P-Wo^i 1 '^)- ( 58 ) 

The asymptotic internal MOND field tp satisfies the Poisson equation for the above asymptotic 
form of p. This determines tp up to a harmonic function (a solution of the Laplace equation). 
Since tp is required to vanish at infinity, only harmonic functions that decrease with r have to be 
considered. Among these we are interested in the leading term, which obviously leaves us with a 
freedom to add to ip a multiple of r _1 . The coefficient of this term is then determined by imposing 
J s Vip-da = 4irGM, which follows by applying Gauss's theorem to eq.(6). Because p oc r~ 3 P2(z/r) 



-17- 



asymptotically, where P2 is the second Legendre polynomial, it follows that the sought after solution 
is of the form r~ 1 P2(z/r). Matching coefficients, we get that 



1 „ MG 3z 2 ^ 

V>* « --vqvq 1 - — 59 

or r A 



solves the above Poisson equation. As can be readily checked, f s Vif)* ■ da = on a sphere. We 
thus have to add to if)* the potential —MG/r to satisfy Gauss's theorem. This finally gives for the 
asymptotic form of tp 

, MG . Oq . - u z 2 

oa + y)! 1 -^^]- (60) 

To compare our results with the analog ones for the nonlinear Poisson formulation, which 
depend on the value of the MOND background acceleration, we have to decide which value of 
go to use when comparing with the results for a given g^. I pick heuristically, the value go = 
z/(|g v |/a )g v . We can then write eq.(58) as 

M Lq 1 3z 2 

p~— A — ttV- 1 -— > 61 

47r/io 1 + Lq r 6 r z 

where po and Lq are, respectively, the values of fi, and of its logarithmic derivative, at rj = |go|/«o 
[fi(x) is related to u{y) as described above]. For the same go, we can write eq.(57) as 

and the asymptotic potential as 

MGl + Lp/2 ( 1 + J«_*\. m 
V ivr 1 + L V 2 + Loi-V 

These can now be compared with the analog quantities for the nonlinear Poisson theory for a 
MOND external field g (Milgrom 2009a): 

M L 1 3z 2 
^-I^Ml + Lo)^ 1 "^' (64) 

where z = (I+Lq) -1 / 2 z, f = (x 2 +y 2 +z 2 ) 1 ^ 2 . We see a similar normalization and radial dependence 
as in p, but a somewhat different angular dependence. For the total "dynamical" mass we have 
there 

1 . ( u ^ 1/2 



M = — arcsin M, (65) 

which, numerically, is not very different from expression (62), differing by at most 6 percent in the 
possible range < Lq < 1. For the asymptotic, internal potential we have 

H r \ 1 + Lq r 2 J 
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An interesting difference between the two theories in the present context is that the symmetry 
axis for the QUMOND theory is the direction of the Newtonian background field, while that in the 
nonlinear Poisson theory is the direction of the MOND background field. The two directions are 
in general different. This shows that the QUMOND theory cannot be equivalent to some form of 
the nonlinear Poisson theory. 16 This difference can also help distinguish observationally between 
the two types of theories. 



5.2. A system dominated everywhere by an external field 

If |Vx| "C I gg' I everywhere in the subsystem, we can write g = go — V^>, with go dictated by 
the mother system, and tp satisfying to lowest order in x 

AV> = "obex* + X, yy + (1 + ^ok, J = 4-irGv p + v v Q x, zz , (67) 
or, in terms of the above, heuristically chosen, go 

AY> = Mo Mx,- + X,w + (1 + £o) -1 X,„]. (68) 

The analog equation for the internal potential in a dominant external field, in the nonlinear Poisson 
formulation is (Milgrom 1986a) 

[^,«« + *P,yy + (1 + Lq)iI>,..\ = —P = PO^X- (69) 

Mo 

In both theories the main effect is to enhance gravity by a factor of l//xo over Newtonian 
gravity. The secondary effect, that of stretching the internal potential in the z direction by a factor 
of about (1 + Lo) 1 / 2 is different in the two theories, as is the direction of the z axis itself. 

5.3. Center-of-mass acceleration of composite systems 

Consider a small body-such as an atom, a gas cloud, or a star-freely falling in a mother system, 
such as a galaxy. Such subsystems are made of constituents that are, sometimes, subject to high 
internal accelerations, hence to high total accelerations. Is it possible then that the center-of-mass 



There is other evidence pointing to the same effect: Consider, as an example, two unequal point masses. As 
explained in Milgrom (1986b), a region where p takes up both signs appears around the critical point where = 
between the two masses. Using the same arguments, we now see that such a region for p appears where V(j> N = 0. 
The critical point for in a nonlinear Poisson theory is never at the Newtonian critical point (unless symmetry 
dictates it). So p cannot be a p for some nonlinear Poisson theory. We shall see more arguments to this effect below. 
If the present theory is equivalent to some nonlinear Poisson formulation with some interpolating function p,, we 
would have to have p(x) = l/u(y), where xp(x) = y in order for the two to coincide for spherical systems. But with 
this identification of v the two theories are definitely not equivalent for aspherical systems. 
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motion of the composite subsystem is subject to the MOND acceleration of the galaxy, as it should 
for MOND phenomenology to work? This was shown to be the case for the nonlinear Poisson 
formulation (Bekenstein & Milgrom 1984), and I now show that it is also the case for the present 
theory. 

Take an arbitrary, bounded mass distribution p, representing the subsystem, placed in a back- 
ground MOND field go, such as a galactic field in the example above. As explained above, provided 
the mass of the subsystem and its extent are much smaller than the corresponding attributes of 
the galaxy, we can describe the MOND field of the system by the solution of the MOND equations 
with a boundary condition of constant acceleration at infinity. 

Use expression (38) to calculate the force on the system as an integral over the surface £. Take 
£ large enough so that on it go strongly dominates the internal acceleration, so we can use the 
asymptotic form of the MOND field g = go — V^>, where ip is given by eq.(60). Substituting this in 
the integrand term — gg^ • da in eq.(38), the term with go gives —go J" g^ • da = 47rGMgo, where 
I used the Newtonian Gauss theorem. The term with VV> combined with all the other terms in 
eq.(38) can be shown to give a vanishing contribution in the limit where S goes to infinity. We are 
then left with F = Mgo, independent of the details of p, giving a center-of-mass acceleration go. 

6. Discussion 

I have presented a new NR formulation of MOND as a modified gravitational potential theory. 
It is derivable from an action and enjoys the standard conservation laws. In a sense, it is an upgrade 
of the pristine formulation of MOND in which the MOND acceleration field is an algebraic function 
of the Newtonian acceleration field. 

The theory is a member in a class of bi-potential theories in which only one-the MOND 
potential-couples to matter directly, whereas the other is an auxiliary potential. The MOND 
departure from Newtonian physics enters not through a modification of the free action of the 
potential (as happens in the original formulation of the Bekenstein Milgrom theory), but through 
the interaction between the two potential fields. In Milgrom (2002b) I described a membrane 
model for MOND in the spirit of the well known membrane model of gravity. In this model MOND 
departure from Newtonian gravity is introduced through modified 'elasticity' of the membrane 
(compared with that corresponding to newtonian gravity). In the same vein, the present class of 
theories gravity can be attributed to the existence of two membranes, one on which matter live 
and with which alone it interacts directly, the other membrane interacting with the first. MOND 
effects are then introduced not as modified elasticity of the membranes-which is normal, albeit 
with the auxiliary membrane having negative elasticity-but through the interaction between the 
membranes. 

These theories are now added to the nonlinear Poisson formulation of MOND propounded a 
quarter century ago by Bekenstein & Milgrom (1984), to which it is similar in many ways. The 
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two theories are shown, in fact, to fork out of the same MOND modification of a Palatini-like 
formulation of Newtonian gravity. 

This addition carries with it the usual benefits of diversity. For example, it brings home the 
possibility that even more formulations await discovering, and encourages us to look for them. 

Also, comparing the predictions of different formulations helps pinpoint results that are not 
generic to the MOND paradigm, not even to its NR formulations. It may also point to new 
directions for constructing relativistic formulations of MOND. 

In itself the new formulation is a complete theory, on a par with the nonlinear Poisson formu- 
lation, but rather easier to work with, involving, as it does, only linear differential equations. 

Because of the similarities between the two formulations, this QUMOND formulation can also 
double as a well motivated approximation for the less wieldy nonlinear Poisson formulation (and 
so also, for the NR limit of theories of the class of TeVeS whose NR limit is a double-potential 
theory of the nonlinear Poisson type). In fact, it is exactly in this role that such a formulation 
was used in Milgrom (2009a)-without recognizing its completeness as a theory-to calculate the 
external-field effect in the inner solar system in the nonlinear Poisson formulation. This yielded a 
closed expression for the desired effect, whereas the application of the nonlinear theory required 
solving the field equation first. Those "approximate" results are now recognized as exact in the 
new theory. 17 



This research was supported by a center of excellence grant from the Israel Science Foundation. 

A. Alternative expressions for the force the torque and the virial 

The force on a subsystem made of the masses within the volume V of surface S, given in 
eq.(33) can also be written as 
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The torque from eq.(39) can be written as 



T = - f d 3 



rpr x V<j> = — 



4ttG 



1 




rA<j) N r x V</>, 



(A2) 



Jv 



17 There were also instances in the past, such as in Milgrom (1986b), and in Milgrom & Sanders (2008) where the 
present formalism was used to approximate the 'phantom' density of the nonlinear Poisson theory; these too are exact 
in the present theory. 
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and similarly for the virial. Integrating by parts several times in a certain order, using the field 
equation (6), and expressing volume integrals of divergences as surface integrals, we get expressions 
(38) (40) (46) for these three quantities. Performing the integration by parts in a different order we 
get other expressions for these quantities, which are also useful 

F = ^{gg" • da + (|q - vgl)da + ^ [(da ■ V)g - V(^ • da)]}, (A3) 

or, writing, as in eq.(8), g = vg N + V x A we get 

F = ~^G Jj^^ ' d<7 + ( f S " ^' )da + (V X A)g ^ ' da + X A) ' da]} - (M) 

T= 4^G Jj~ T X g(gJV ' dCT) " f Qr X dCT + ^ X ' dCT)+ 

+<^[V • (isg N )(r x da) + g x da - (r x V)(g • da)]}. (A5) 

AttGV = J d\\W N - \a\Q\ + J r • g(g" ■ da) + ± J [a 2 Q - 2vg 2 N ]r • da+ 

+ J <f> N [(g - i/g") ■ da + (r • V)(g - i/g") • da]. (A6) 
In all the above expressions V in the last term does not act on the surface element da. 
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